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Upper limit map of a background of gravitational waves 
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We searched for an anisotropic background of gravitational waves using data from the LIGO S4 
science run and a method that is optimized for point sources. This is appropriate if, for example, 
the gravitational wave background is dominated by a small number of distinct astrophysical sources. 
No signal was seen. Upper limit maps were produced assuming two different power laws for the 
source strain power spectrum. For an /~ 3 power law and using the 50 Hz to 1.8kHz band the 
upper limits on the source strain power spectrum vary between 1.2 x 10 _4S Hz~ 1 (100 Hz//) 3 and 
1.2 x 10 _47 Hz _1 (100 Hz//) 3 , depending on the position in the sky. Similarly, in the case of constant 
strain power spectrum, the upper limits vary between 8.5 x 10 _49 Hz _1 and 6.1 x 10 _48 Hz _1 . As a 
side product a limit on an isotropic background of gravitational waves was also obtained. All limits 
are at the 90% confidence level. Finally, as an application, we focused on the direction of Sco-Xl, 
the closest low-mass X-ray binary. We compare the upper limit on strain amplitude obtained by 
this method to expectations based on the X-ray luminosity of Sco-Xl. 

PACS numbers: 04.80.Nn, 04.30.Db, 07.05.Kf, 02.50.Ey, 02.50.Fz, 95.55.Ym, 98.70.Vc 



I. INTRODUCTION 

A stochastic background of gravitational waves can 
be non-isotropic if, for example, the dominant source of 
stochastic gravitational waves comes from an ensemble of 
astrophysical sources (e.g. d,[IH), and if this ensemble 
is dominated by its strongest members. So far the LIGO 
Scientific Collaboration has analyzed the data from the 
first science runs for a stochastic background of gravita- 
tional waves [EH i], assuming that this background is 
isotropic. If astrophysical sources indeed dominate this 
background, one should look for anisotropies. 

A method that is optimized for extreme anisotropies, 
namely point sources of stochastic gravitational radia- 
tion, was presented in Q- It is based on the cross- 
correlation of the data streams from two spatially sepa- 
rated gravitational wave interferometers, and is referred 
to as radiometer analysis. We have analyzed the data of 
the 4th LIGO science run using this method. 

Section |TT] is a short description of the radiometer anal- 
ysis method. The peculiarities of the S4 science run are 
summarized in section UITl and we discuss the results in 
section [TV] 



II. METHOD DESCRIPTION 

A stochastic background of gravitational waves can 
be distinguished from other sources of detector noise 
by cross-correlating two independent detectors. Thus 
we cross-correlate the data streams from a pair of de- 
tectors with a cross-correlation kernel Q, chosen to be 
optimal for a source which is specified by an assumed 
strain power spectrum H (/) and angular power distribu- 
tion P(f2). Specifically, with si(/) and S2(/) representing 
the Fourier transforms of the strain outputs of two detec- 
tors, this cross-correlation is computed in the frequency 
domain segment by segment as: 



Yt = 



dfSl(f)Qt(f)Uf). 



(1) 



In contrast to the isotropic analysis the optimal filter Q t 
is now sidereal time dependent. It has the general form: 



HI) = At 



fs* dCl 1(lit (f)P(Cl)H(f) 



(2) 



where At is a normalization factor, P± and P2 are the 
strain noise power spectra of the two detectors, H is the 
strain power spectrum of the stochastic signal we search 
for, and the factor 7^ t takes into account the sidereal 
time dependent time delay due to the detector separation 
and the directionality of the acceptance of the detector 
pair. Assuming that the source is unpolarized, 7^ t is 
given by: 



7ft. 



1 



t (/) = oE 



(3) 



where /S.x = x 2 — £1 is the detector separation vector, 
Cl is the unit vector specifying the sky position and 



e ab 



yayb 



(4) 



is the response of detector i to a zero frequency, unit am- 
plitude, A — + or x polarized gravitational wave. e^ b (fl) 
is the spin-two polarization tensor for polarization A and 
Xf and Y°" are unit vectors pointing in the directions of 
the detector arms (see [f| for details). The sidereal time 
dependence enters through the rotation of the earth, af- 
fecting Xf, Y? and Ax. 

The optimal filter Q t is derived assuming that the in- 
trinsic detector noise is Gaussian and stationary over 
the measurement time, uncorrelated between detectors, 
and uncorrelated with and much greater in power than 
the stochastic gravitational wave signal. Under these 
assumptions the expected variance, o~\ , of the cross- 
correlation is dominated by the noise in the individ- 
ual detectors, whereas the expected value of the cross- 
correlation Yt depends on the stochastic background 
power spectrum: 



(Y t 2 ) - (Y t f 



T 



{Qt, Qt) 



(5) 
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(Yt) = T 



Ss> dn lfltt P((l)H s 

P1P2 



(0) 



Here the scalar product (-,-) is denned as (A, B) = 
/ M oo A*(/)B(/)F 1 (/)P 2 (/)d/ and T is the duration of 
the measurement. 

Equation [5] defines the optimal filter Q t for any arbi- 
trary choice of P(f2). To optimize the method for finite, 
but unresolved astrophysical sources one should use a 
P(f2) that covers only a localized patch in the sky. But 
the angular resolution is diffraction limited with the de- 
tector separation as baseline and the frequency content 
weighted by i? 2 P 1 _1 P 2 ~ 1 . For a constant H(f) this corre- 
sponds to a resolution of several tens of square degrees, so 
that astrophysical sources will not be spatially resolved. 
Thus we chose to optimize the method for true point 
sources, i.e. P(O) = 8 2 (Cl,tt'), which also allows for a 
more efficient implementation (see Q). 

We define the strain power spectrum H(f) of a point 
source as one-sided (positive frequencies only) and in- 
cluding the power in both polarizations. Thus H (/) is re- 
lated to the gravitational luminosity Lqw and the grav- 
itational energy flux F(f) through 



L 



GW 



fx, 



/min 



H(f)fdf, (7) 



with c the light speed and G Newton's constant. We look 
for strain power spectra H(f) in the form of a power law 
with exponent (3. The amplitude at the pivot point of 
100 Hz is described by Hp, i.e. 



H(f) = Hp 



f 



100 Hz 



(8) 



With this definition we can choose the normalization of 
the optimal filter Q t such that equation [5] reduces to 



(Yt) = H t 



0- 



(9) 



The data set from a given interferometer pair is divided 
into equal-length intervals, and the cross-correlation Y t 
and theoretical oy t are calculated for each interval, yield- 
ing a set {Y t ,<TY t } of such values for each sky direction 
Q, with t the mid-segment sidereal time. The optimal 
filter Q t is kept constant and equal to its mid-segment 
value for the whole segment. The remaining error due to 
this discretization is of second order in (T seg /1 day) and 
is given by: 



^err (Ps 



SCg 




(10) 



with / the typical frequency and d the detector separa- 
tion. At the same time the interval length can be cho- 
sen such that the detector noise is relatively stationary 



over one interval. We use an interval length of 60 sec, 
which guarantees that the relative error Y erT (T seg )/Y is 
less than 1%. The cross-correlation values are combined 
to produce a final cross-correlation estimator, Y op f, that 
maximizes the signal-to-noise ratio, and has variance 



'opt- 
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(11) 



In practice the intervals are overlapping by 50% to avoid 
the effective loss of half the data due to the required 
windowing (Hanning) . Thus equation [TT] was modified 
slightly to take the correlation of neighboring segments 
into account. 

The data was downsampled to 4096 Hz and high-pass 
filtered with a sixth order Butterworth filter with a cut- 
off frequency at 40 Hz. Frequencies between 50 Hz and 
1800 Hz were used for the analysis and the frequency 
bin width was 0.25 Hz. Frequency bins around multi- 
ples of 60 Hz up to the tenth harmonic were removed, 
along with bins near a set of nearly monochromatic in- 
jected signals used to simulate pulsars. These artificial 
pulsars proved useful in a separate end-to-end check of 
this analysis pipeline, which successfully recovered the 
sky locations, frequencies and strengths of three such 
pulsars listed in TABLE fl] The resulting map for one 
of these pulsars is shown in FIG. [5] 



III. THE LIGO S4 SCIENCE RUN 

The LIGO S4 science run consisted of one month of co- 
incident data taking with all three LIGO interferometers 
(22 Feb 2005 noon to 23 Mar 2005 midnight CST). Dur- 
ing that time all three interferometers where roughly a 
factor of 2 in amplitude away from design sensitivity over 
almost the whole frequency band. Also, the Livingston 
interferometer was equipped with a Hydraulic External 
Pre-Isolation (HEPI) system, allowing it to stay locked 
during day time. This made S4 the first LIGO science 
run with all-day coverage at both sites. A more detailed 
description of the LIGO interferometers is given in [5j . 

Since the radiometer analysis requires two spatially 
separated sites we used only data from the two 4 km 
interferometers (HI in Hanford and LI in Livingston). 
For these two interferometers about 20 days of coinci- 
dent data was collected, corresponding to a duty factor 
of 69%. 

The large spatial separation also reduces environmen- 
tal correlations between the two sites. Nevertheless we 
still found a comb of 1 Hz harmonics that was coherent 
between HI and LI. This correlation was found to be 
at least in part due to an exactly 1-sec periodic signal 
in both interferometers (FIG. [T|), which was caused by 
cross-talk from the GPS_RAMP signal. The GPS.RAMP 
signal consists of a 10 msec saw-tooth signal that starts 
at every full second, lasts for 1 msec and is synchronized 
with the GPS receivers, (see FIG. [IJ. This ramp was 
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Injected pulsars 

Quantity 


Pulsar #3 


Pulsar #4 


Pulsar #8 


Frequency during S4 run 


108.86 Hz 


1402.20 Hz 


193.94 Hz 


Noise level (er) 


1.89 x 1CT 47 


6.04 x 10 -46 


1.73 x 10 -47 


Injected Hdf (corrected for polarization) 


1.74 x 1(T 46 


4.28 x 10~ 44 


1.54 x 10~ 46 


Recovered Hdf on source 


1.74 x 1CT 46 


4.05 x 10 -44 


1.79 x 10~ 46 


Qi cm q ci- ti m cp T"fti"i <~\ ( KlSJ W 1 


9.2 


67.1 


10.3 


Injected position 


llh 53m 29.4s 


18h 39m 57.0s 


23h 25m 33.5s 




-33° 26' 11.8" 


-12° 27' 59.8" 


-33° 25' 6.7" 


Recovered position (max SNR) 


12h 12m 


18h 40m 


23h 16m 




-37° 


-13° 


-32° 



TABLE I: Injected pulsars: The table shows the level at which the three strongest injected pulsars were recovered. Hdf 
denotes the RMS strain power over the 0.5 Hz band that was used. The reported values of for the injected Hdf include 
corrections that account for the difference between the polarized pulsar injection and an unpolarized source that is expected 
by the analysis. The underestimate of Pulsar #4 is due to a known bias of the analysis method in the case of a signal strong 
enough to affect the power spectrum estimation. 



used as an off-line monitor of the Analog-to-Digital Con- 
verter (ADC) card timing and thus was connected to the 
same ADC card that was used for the gravitational wave 
channel, which resulted in a non-zero cross-talk to the 
gravitational wave channel. 

To reduce the contamination from this signal a tran- 
sient template was subtracted in the time domain. This 
has the advantage that effectively only a very narrow 
band (1/runtime ~ 1 x 1CP 6 Hz) is removed around 
each 1 Hz harmonic, while the rest of the analysis is 
unaffected. The waveform for subtraction from the raw 
(uncalibrated) data was recovered by averaging the data 
from the whole run in order to produce a typical second. 
Additionally, since this typical second only showed signif- 
icant features in the first 80 msec, the transient subtrac- 
tion template was set to zero (with a smooth transition) 
after 120 msec. This subtraction was done for only HI 
since adding repetitive data to both detectors can intro- 
duce artificial correlation. It eliminated the observed cor- 
relation. However, due to an automatically adjusted gain 
between the ADC card and the gravitational wave chan- 
nel, the amplitude of the transient waveform is affected 
by a residual systematic error. Its effect on the cross- 
correlation result was estimated by comparing maps with 
the subtraction done on either HI or LI. The systematic 
error is mostly concentrated around the north and south 
poles, with a maximum of about 50% of the statistical 
error at the south pole. In the declination range of —75° 
to +75° the error is less than 10% of the statistical er- 
ror. For upper limit calculations this systematic error is 
added in quadrature to the statistical error. After the S4 
run the GPS_RAMP signal was replaced with a two-tone 
signal at 900 Hz and 901 Hz. The beat between the two 
is now used to monitor the timing. 

One post-processing cut was required to deal with de- 
tector non-stationarity. To avoid a bias in the cross- 
correlation statistics the segment before and the segment 
after the one being analyzed are used for the power spec- 
tral density (PSD) estimate [T(ij |. Therefore the analysis 



becomes vulnerable to large, short transients that happen 
in one instrument in the middle segment - such transients 
cause a significant underestimate of the PSD and thus of 
the theoretical standard deviation for this segment. This 
leads to a contamination of the final estimate. 

To eliminate this problem the standard deviation a 
is estimated for both the middle segment and the two 
adjacent segments. The two estimates are then required 
to agree within 20%: 

— < (Jmiddlc < 1.2. (12) 

1 • 2 ^"adj acent 

The analysis is fairly insensitive to the threshold - the 
only significant contamination comes from very large out- 
liers that are cut by any reasonable threshold [8]. The 
chosen threshold of 20% eliminates less than 6 % of the 
data. 



IV. RESULTS FROM THE S4 RUN 
A. Broadband results 

In this analysis we searched for an H{f) following a 
power law with two different exponents (3: 

. P = -3: H(J)=H- 3 (^Y*) 3 . 
This emphasizes low frequencies and is useful when 
interpreting the result in a cosmological framework, 
since it corresponds to a scale-invariant primordial 
perturbation spectrum, i.e. the GW energy per 
logarithmic frequency interval is constant. 

• = 0: H(f) = Hq (constant strain power). 
This emphasizes the frequencies for which the in- 
terferometer strain sensitivity is highest. 

The results are reported as point estimate and cor- 
responding standard deviation er^ for each pixel Q. The 
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H1:DARM_ERR (1437280 averages) 




100 
msec 

H1:DARM_ERR (1437280 averages) 




100 
msec 

L1:DARM_ERR (1447904 averages) 




FIG. 1: Periodic timing transient in the gravitational 
wave channel (D ARMJ3RR) , calibrated in //Volt at the ADC 
(Pentek card) for HI (top two) and LI (bottom two) shown 
with a span of 200 msec and 14 msec in black. The x-axis is 
the offset from a full GPS second. About 1.4 million seconds 
of DARMJ3RR data was averaged to get this trace. Also 
shown in gray is the GPS_RAMP signal that was used as a 
timing monitor. It was identified as a cause of the periodic 
timing transient in DARM_ERR. The HI trace shows an ad- 
ditional feature at 6 msec. 



point estimate Ya must be interpreted as best fit ampli- 
tude Hp for the pixel Cl (equation [9|) . 

Also we should note that the resulting maps have an 
intrinsic spatial correlation, which is described by the 
point spread function 



A(Cl,Cl') = 



(13) 



It describes the spatial correlation in the following sense: 
if cither Y^, — Y due to random fluctuations, or if there 
is a true source of strength Y at Cl', then the expec- 
tation value at Cl is (Y^) = A(Cl,Cl')Y. The shape of 

A(fi, Cl') depends strongly on the declination. FIG. 0] 
shows A(Cl, Cl') for different source declinations and both 
the [3 = —3 and /3 = case, assuming continuous day 
coverage. 



Scale-invariant case, f3 = —3 



S4, beta=-3 Histogram of SNR (40 bins) 




Data 

Ideal gaussian (sigma=1 mean=0) 

Max Likelihood: sigma=0.91836 mean=0.11816 

1 -sigma error for 100 indep. points 



FIG. 2: S4 Result: Histogram of the signal-to-noise ratio 
(SNR) for /3 = —3. The gray curve is a maximum likeli- 
hood Gaussian fit to the data. The black solid line is an 
ideal Gaussian, the two dash-dotted black lines indicate the 
expected one sigma variations around this ideal Gaussian for 
100 independent points (N e s = 100). 

A histogram of the SNR = £ i s plotted in FIG.0 The 
data points were weighted with the corresponding sky 
area in square degrees. Because neighboring points are 
correlated, the effective number of independent points 
N c g is reduced. Therefore the histogram can exhibit 
statistical fluctuations that are significantly larger than 
those naively expected from simply counting the number 
of pixels in the map, while still being consistent with (cor- 
related) Gaussian noise. Indeed the histogram in FIG. [2] 
features a slight bump around SNR=2, but is still con- 
sistent with N c ff — 100 - the red dash-dotted lines indi- 
cate the one sigma band around the red ideal Gaussian 
for Nes = 100. Additionally the SNR distribution also 
passes a Kolmogorov-Smirnov test for N a g = 100 at the 
90% significance level. 

The number of independent points N c s , which in effect 
describes the diffraction limit of the LIGO detector pair, 
was estimated by 2 heuristic methods: 

• Spherical harmonics decomposition of the SNR 
map. The resulting power vs I graph shows struc- 
ture up to roughly I = 9 and falls off steeply above 
that - the I = 9 point corresponds to one twenti- 
eth of the maximal power. The effective number of 
independent points then is N e g ~ (I + l) 2 = 100. 

• FWHM area of a strong injected source, which is 
latitude dependent but of the order of 800 square 
degrees. To fill the sky we need about N c g w 50 of 
those patches. We used the higher estimate N e g = 
100 for this discussion. 
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FIG. [2] suggests that the data is consistent with no sig- 
nal. Thus we calculated a Bayesian 90% upper limit for 
each sky direction. The prior was assumed to be flat be- 
tween zero and an upper cut-off set to 5 x 10 _45 Hz _1 at 
100 Hz, the approximate limit that can be set from just 
operating a single LIGO interferometer at the S4 sensi- 
tivity. Note, however, that this cut-off is so high that 
the upper limit is completely insensitive to it. Addition- 
ally we marginalized over the calibration uncertainty of 
8% for HI and 5% for LI using a Gaussian probability 
distribution. The resulting upper limit map is shown 
in FIG. [HI The upper limits on the strain power spec- 
trum H(f) vary between 1.2 x lO" 48 ^" 1 (100 Hz//) 3 
and 1.2 x lO-^Hz" 1 (100 Hz//) 3 , depending on the po- 
sition in the sky. These strain limits correspond to 
limits on the gravitational wave energy flux F(f) vary- 
ing between 3.8 x 10~ 6 erg cm" 2 Hz~ 4 (100 Hz//) and 
3.8 x 10" 5 erg cm" 2 Hz" 1 (100 Hz//). 



2. Constant strain power, (3 — 

Similarly, FIG. [3] shows a histogram of the SNR = — 
for the constant strain power case. Structure in the 
spherical harmonics power spectrum goes up to I = 19, 
thus N e s was estimated to be N cfi » (I + l) 2 = 400. Al- 
ternatively the FWHM area of a strong injection covers 
about 100° 2 which also leads to N e g w 400. The dash- 
dotted red lines in the histogram (FIG. [3]) correspond to 
the expected 1 — a deviations from the ideal Gaussian for 
N c g = 400. The histogram is thus consistent with (cor- 
related) Gaussian noise, indicating that there is no signal 
present. The SNR distribution also passes a Kolmogorov- 
Smirnov test for N e g = 400 at the 90% significance level. 

Again we calculated a Bayesian 90% upper limit for 
each sky direction, including the marginalization over the 
calibration uncertainty. The prior was again assumed to 
be flat between and an upper cut-off of 5 x 10 _45 Hz _1 
at 100 Hz. The resulting upper limit map is shown in 
FIG. [7J The upper limits on the strain power spectrum 
H(f) vary between 8.5 x lO^Hz -1 and 6.1 x 10- 48 Hz _1 
depending on the position in the sky. This corresponds to 
limits on the gravitational wave energy flux F(f) varying 
between 2.7 x 10~ 6 erg cm~ 2 Hz -1 (//100 Hz) 2 and 1.9 x 
10" 5 erg cm" 2 Hz" 1 (//100 Hz) 2 . 



3. Interpretation 

The maps presented in FIGS. [5] and [7] represent the 
first directional upper limits on a stochastic gravitational 
wave background ever obtained. They are consistent with 
no gravitational wave background being present. This 
search is optimized for well localized, broadband sources 
of gravitational waves. As such it is best suited for un- 
expected, poorly modeled sources. 



S4, beta=0 Histogram of SNR (40 bins) 




J Data 

- Ideal gaussian (sigma=1 mean=0) 
Max Likelihood: sigma=0.99738 mean=-0. 025485 
1 -sigma error for 400 indep. points 



FIG. 3: S4 Result: Histogram of the signal-to-noise ratio 
(SNR) for P — 0. The gray curve is a maximum likelihood 
Gaussian fit to the data. The black solid line is an ideal 
Gaussian, the two dash-dotted black lines indicate the ex- 
pected one sigma variations around this ideal Gaussian for 
400 independent points (JV e g = 400). 



In order to compare the result to what could be ex- 
pected from known sources we also search for the gravita- 
tional radiation from low-mass X-ray binaries (LMXBs). 
They are accretion-driven spinning neutron stars, i.e. 
narrow-band sources and thus not ideal for this broad- 
band search. However they have the advantage that we 
can predict the gravitational luminosity based on the 
known X-ray flux. If gravitational radiation provides the 
torque balance for LMXBs, then there is a simple relation 
between the grav itational luminosity Lqw and X-ray lu- 
minosity Lx [13| : 



Lew 



f 



spin 



/Kepi 



Lx- 



(14) 



Here /Kepler is final orbital frequency of the accreting 
matter, about 2 kHz for a neutron star, and / sp i n is the 
spin frequency. 

As an example we estimate the gravitational luminos- 
ity of all LMXBs within the Virgo galaxy cluster. Their 
integrated X-ray luminosity is about 10 -9 erg/sec/cm 2 
(3000 galaxies at 15 Mpc, 10 40 erg/sec/galaxy from 
LMXBs). For simplicity we assume that the ensemble 
produces a flat strain power spectrum H(f) over a band- 
width A/. Then the strength of this strain power spec- 
trum is about 



H(f) 



2G 



TTC"^ j "Kepler j center^ f 



-L 



x 



10" 55 Hz~ 



100 Hz 

/center 



100 Hz 



(15) 



Here /center is the typical frequency of the A/ wide band 
of interest. This is quite a bit weaker than the upper limit 
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set in this paper, which is mostly due to the fact that the 
intrinsically narrow-band sources are diluted over a broad 
frequency band. 
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FIG. 4: Point spread function A(Cl, Q') of the radiome- 
ter for P = —3 (top two figures) and for f3 = (bottom two 
figures). Plotted is the relative expected signal strength as- 
suming a source at right ascension 12 h and declinations 20° 
and 60°. Uniform day coverage was assumed, so the resulting 
shapes are independent of right ascension. An Aitoff projec- 
tion was used to plot the whole sky. 




Right ascension [hours] 

FIG. 5: Injected pulsar #3: The analysis was run using the 
108.625 Hz — 109.125 Hz frequency band. The artificial signal 
of Pulsar #3 at 108.86 Hz stands out with a signal-to-noise 
ratio of 9.2. The circle marks the position of the simulated 
pulsar. 



S4, beta=-3 90% confidence upper limit 
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FIG. 6: S4 Result: Map of the 90% confidence level 
Bayesian upper limit on Hp for /3 = —3. The upper 
limit varies between 1.2 x 10" 48 Hz _1 (100 Hz//) 3 and 1.2 x 
10 _47 Hz~ 1 (100 Hz//) 3 , depending on the position in the sky. 
All fluctuations are consistent with the expected noise. 



S4, beta=0 90% confidence upper limit 
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B. Limits on isotropic background 



FIG. 7: S4 Result: Map of the 90% confidence level 
Bayesian upper limit on Hp for /3 = O.The upper limit varies 
between 8.5 x 10" 49 Hz _1 and 6.1 x 10" 48 Hz _1 depending on 
the position in the sky. 



It is possible to recover the estimate for an isotropic 
background as an integral over the map (see 0]). The cor- 
responding theoretical standard deviation would require 
a double integral with essentially the point spread func- 
tion as integrand. In practice it is simpler to calculate 
this theoretical standard deviation directly by using the 



overlap reduction function for an isotropic background 
(see pi]). From that the 90% Bayesian upper limit can 
be calculated, which is additionally marginalized over the 
calibration uncertainty. In the j3 = —3 case the 90% up- 
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per limit we can set on h^ 2 ^igw(f) is 1-20 x 10~ 4 . The 
dimensionless quantity fi gw (/) is the GW energy density 
per unit logarithmic frequency, divided by the critical 
energy density p c to close the universe, and /172 is the 
Hubble constant in units of 72 km sec" ^Mpc -1 . TablcHT] 
summarizes the results for all choices of (3. 



1. Interpretation 

In [l[ we published an upper limit of h^ 2 ^GW < 
6.5 x 10~ 5 on an isotropic gravitational wave background 
using S4 data. That analysis is mathematically identi- 
cal to inferring the point estimate as an integral over the 
map Q, but the mitigation of the timing transient and 
the data quality cuts were sufficiently different to affect 
the point estimate. While both results are consistent 
within the error bar of the measurement, this difference 
results in a slightly higher upper limit. Both results are 
significantly better than the previously published LIGO 
S3 result. 



C. Narrow-band results targeted on Sco-Xl 

As an application we again focus on low-mass X-ray 
binaries (LMXBs). The gravitational wave flux from all 
LMXBs is expected to be dominated by the closest one, 
Sco-Xl, simply because Sco-Xl dominates the X-ray flux 
from all LMXBs, and X-ray luminosity Lx is related to 
the gravitational luminosity Lqw through equation 1141 
Unfortunately the spin frequency of Sco-Xl is not known. 
We thus want to set an upper limit for each frequency bin 
on the RMS strain coming from the direction of Sco-Xl 
(RA: 16h 19m 55.0850s; Dec: -15° 38' 24.9"). 

The binary orbital velocity of Sco-Xl is about 40 ± 
5 km/sec (see [l2j). This induces a maximal frequency 
shift of Af GW = 2.7 x 10~ 4 x few- We chose a bin 
width of df = 0.25 Hz, which is broader than maximal 
frequency shift Afcw for all frequencies below 926 Hz 
and is the same bin width that was used for the broad- 
band analysis. Above 926 Hz the signal is guaranteed to 
spread over multiple bins. 

To avoid contamination from the hardware-injected 
pulsars, the 2 frequency bins closest to a pulsar frequency 
were excluded. Multiples of 60 Hz were also excluded. 
The lowest frequency bin was at 50 Hz, the highest one 
at 1799.75 Hz. FIG. [5] shows a histogram of the remain- 
ing 6965 0.25 Hz wide frequency bins. It is consistent 
with a Gaussian distribution (Kolmogorov-Smirnov test 
with N — 6965 at the 90% significance level). 

A 90% Bayesian upper limit for each frequency bin 
was calculated based on the point estimate and standard 
deviation, including a marginalization over the calibra- 
tion uncertainty. Figure [5] is a plot of this 90% limit (red 
trace) . Above about 200 Hz (shot noise regime above the 
cavity pole) the typical upper limit rises linearly with fre- 
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FIG. 8: S4 Result for Sco-Xl: Histogram of the signal- 
to-noise ratio calculated for each 0.25 Hz wide frequency bin. 
There are no outliers. 



Sco-X1, ra= 16.332 h , decl= -15.6402 deg 
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FIG. 9: S4 Result for Sco-Xl: The 90% confidence 
Bayesian upper limit as a function of frequency - marginalized 
over the calibration uncertainty. The standard deviation (one 
sigma error bar) is shown in blue. 



quency and is given by 

hS « 3.4 x 10- 24 (2^) / > 200 Hz. (16) 
The standard deviation is also shown in blue. 

1. Interpretation 

In principle, the radiometer analysis is not an optimal 
method to search for a presumably periodic source like 
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S4 isotropic upper limit 

Quantity 


«gw(/) = 


const 




Sgw (/) = const 


point estimate Y 
standard deviation a 


1.02 x 10" 
6.97 x 10" 


-A7-T-T - 

47 Hz 
- 48 Hz- 


\ 100Hz ) „ 

linnHz/ 


-7.12 x W^Uz^ 1 
7.22 x lO^^Hz^ 1 


90% Bayesian UL on S gw (/) 
90% Bayesian UL on h^^gwif) 


1.99 x 10" 
1.20 x 10" 


- 47 Hz" 

-4 


V 100Hz ) 


8.49 x lO-^Hz" 1 

5-13 x io- 5 



TABLE II: S4 isotropic result for the f2 gw (/) = const, (f3 = —3) and the 5 , gw (/) = const, (/3 = 0) case. The first two lines 
show the point estimate and standard deviation that are used to calculate the 90% Bayesian upper limits. The upper limits 
are also marginalized over the calibration uncertainty. These results agree with the ones published in [l|] within the error bar 
of the measurement. 



Sco-Xl. Nevertheless it can set a competitive upper limit 
with a minimal set of assumptions on the source and 
significantly less computational resources. Indeed LIGO 
published a 95% upper limit on gravitational radiation 
amplitude from Sco-Xl of 1.7xl0~ 22 to 1.3xl0~ 21 across 
the 464 - 484 Hz and 604 - 624 Hz frequency bands 0], 
using data from S2, which had a noise amplitude about 
4.5 times higher around 500 Hz in each instrument. The 
analysis was computationally limited to using 6 hours of 
data and two 20 Hz frequency bands. However the strain 
amplitude sensitivity scales as J 1-1 / 4 Q j while a coherent 
method scales as T~ x /' 2 . 

The upper limit (eq. [Tfj)) can directly be compared to 
the expected strain based on the X-ray luminosity: 



, (90%) 
RMS 

RMS 



100 



200 Hz 



/> 200 Hz. (17) 



Here / is the gravitational wave frequency, i.e. twice 
the (unknown) spin frequency of Sco-Xl. This is close 
enough that, if the model described in and thus 

equation Q4] are indeed correct, Sco-Xl ought to be de- 
tectable with this method and the next generation of 
gravitational wave detectors operated in a narrow-band 
configuration (AdvLIGO [3). 

For a discussion of the 
expected signal from Sco-Xl see also 



V. CONCLUSION 

We produced the first upper limit maps for a stochas- 
tic gravitational wave background by applying a method 
that is described in [7J to the data from the LIGO S4 sci- 
ence run. No signal was seen and upper limits were set for 
two different choices for the strain power spectrum H(f). 
In the case of H(f) cx /~ 3 the upper limits for a point 
source vary between 1.2 x lO-^Hz" 1 (100 Hz//) 3 and 
1.2 x 10 _47 Hz~ 1 (100 Hz//) 3 , depending on the position 
in the sky (see FIG. [B]). Similarly, in the case of constant 
H(f) the upper limits vary between 8.5 x 10 _49 Hz _1 and 
6.1 x 10 _48 Hz _1 (see FIG. [7|). As a side product limits 



on an isotropic background of gravitational waves were 
also obtained, see TABLE [H] 

In an additional application, narrow-band upper lim- 
its were set on the gravitational radiation coming from 
the closest low- mass X-ray binary, Sco-Xl (see FIG. [9]). 
In the shot noise limited frequency band (above about 
200 Hz) the limits on the strain in each 0.25 Hz wide 
frequency bin follow roughly 



hS « 3.4 x IO- 4 



( f 



V200 Hz 



/> 200 Hz, (18) 



where / is the gravitational wave frequency (twice the 
spin frequency). 
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